Scalable  Coupling  of  Multiscale  AEH  and  PARADYN 
Analyses  for  Impact  Modeling 


by  Rama  R.  Valisetty,  Peter  W.  Chung,  and  Raju  R.  Namburu 


ARL-TR-3512 


June  2005 


Approved  for  public  release;  distribution  is  unlimited. 


NOTICES 

Disclaimers 

The  findings  in  this  report  are  not  to  be  construed  as  an  official  Department  of  the  Army  position  unless 
so  designated  by  other  authorized  documents. 

Citation  of  manufacturer’s  or  trade  names  does  not  constitute  an  official  endorsement  or  approval  of  the 
use  thereof. 


Destroy  this  report  when  it  is  no  longer  needed.  Do  not  return  it  to  the  originator. 


Army  Research  Laboratory 

Aberdeen  Proving  Ground,  MD  21005-5067 


ARL-TR-3512 


June  2005 


Scalable  Coupling  of  Multiscale  AEH  and  PARADYN 
Analyses  for  Impact  Modeling 

Rama  R.  Valisetty,  Peter  W.  Chung,  and  Raju  R.  Namburu 
Computational  and  Information  Sciences  Directorate,  ARL 


Approved  for  public  release;  distribution  is  unlimited. 


Report  Documentation  Page 


Form  Approved 
OMB  No.  0704-0188 


Public  reporting  burden  for  this  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data  sources,  gathering 
and  maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  information.  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of 
information,  including  suggestions  for  reducing  the  burden,  to  Department  of  Defense,  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports  (0704-0188), 
1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington,  VA  22202-4302.  Respondents  should  be  aware  that  notwithstanding  any  other  provision  of  law,  no  person  shall  be  subject  to  any  penalty 
for  failing  to  comply  with  a  collection  of  information  if  it  does  not  display  a  currently  valid  OMB  control  number. 

PLEASE  DO  NOT  RETURN  YOUR  FORM  TO  THE  ABOVE  ADDRESS. 


1.  REPORT  DATE  (DD-MM-YYYY) 

June  2005 


4.  TITLE  AND  SUBTITLE 


2.  REPORT  TYPE 


Scalable  Coupling  of  Multiscale  AEH  and  PARADYN  Analyses  for  Impact 
Modeling 


3.  DATES  COVERED  (From  -  To) 

2003-2004 


5a.  CONTRACT  NUMBER 


5b.  GRANT  NUMBER 


5c.  PROGRAM  ELEMENT  NUMBER 


6.  AUTHOR(S) 

Rama  R.  Valisetty,  Peter  W.  Chung,  and  Raju  R.  Namburu 


5d.  PROJECT  NUMBER 

5UH7FL 


5e.  TASK  NUMBER 


5f.  WORK  UNIT  NUMBER 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES) 

U.S.  Army  Research  Laboratory 

ATTN:  AMSRD-ARL-CI-HC 

Aberdeen  Proving  Ground,  MD  21005-5067 


9.  SPONSORING/MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 


8.  PERFORMING  ORGANIZATION 
REPORT  NUMBER 

ARL-TR-3512 


10.  SPONSOR/MONITOR'S  ACRONYM(S) 


12.  DISTRIBUTION/AVAILABILITY  STATEMENT 

Approved  for  public  release;  distribution  is  unlimited. 


11.  SPONSOR/MONITOR'S  REPORT 
NUMBER(S) 


14.  ABSTRACT 

This  report  describes  scalable  coupling  of  two  stand-alone  computer  codes  for  a  multiscale  impact  analysis  of  composites.  An  asymptotic  expansion 
homogenization  (AEH)-based  microstructural  model  available  for  modeling  microstructural  aspects  of  modem  armor  materials  is  coupled  with 
PARADYN,  a  parallel  explicit  Lagrangian  finite-element  code.  The  first  code  enables  modeling  of  a  material  microstructure  and  provides  material 
response  in  terms  of  global  structural  response  at  a  material  integration  point.  Microstructural  codes  such  as  this  one  are  typically  used  in  stand-alone 
form  and  applied  in  simple  loading  situations.  The  coupling  provided  here  enables  a  micro/macro  type  multiscale  analysis  under  generalized 
three-dimensional  loading  conditions.  Three  sets  of  results  are  presented  to  demonstrate:  (1)  the  verification  of  the  AEH-PARADYN  model 
coupling  to  PARADYN,  (2)  the  scalability  of  the  coupled  model,  and  (3)  an  application  for  modeling  the  impact  response  of  armor  materials. 

In  the  present  work,  a  consistent  AEH  numerical  fonnulation  was  selected  that  was  earlier  shown  to  be  scalable  on  different  computing  architectures. 
In  conjunction  with  a  second-order  accurate  velocity-based  explicit  time  integration  method,  the  formulation  is  coupled  within  the  message-passing 
interface  scheme  of  PARADYN.  PARADYN  is  a  scalable  version  of  the  Lawrence  Livermore  National  Laboratory's  serial  DYNA3D  explicit 
Lagrangian  finite-element  code  used  for  obtaining  large  deformation,  elastic/plastic  response  of  structures. 


15.  SUBJECT  TERMS 

AEH,  PARADYN,  HPC,  scalable  computing,  impact  modeling 


16.  SECURITY  CLASSIFICATION  OF: 

17.  LIMITATION 

OF  ABSTRACT 

18.  NUMBER 

OF  PAGES 

19a.  NAME  OF  RESPONSIBLE  PERSON 

Rama  R.  Valisetty 

a.  REPORT 

UNCLASSIFIED 

b.  ABSTRACT 

UNCLASSIFIED 

c.  THIS  PAGE 

UNCLASSIFIED 

UL 

44 

19b.  TELEPHONE  NUMBER  ( Include  area  code) 

410-278-6638 

Standard  Form  298  (Rev.  8/98) 

Prescribed  by  ANSI  Std.  Z39.18 


Contents 


List  of  Figures  iv 

List  of  Tables  v 

1.  Introduction  1 

2.  AEH  Development  3 

3.  Coupling  the  Microstructural  AEH  Model  in  PARADYN  4 

3 . 1  Guidelines  for  U sing  the  Model . 4 

3.2  Details  of  the  Data  Exchange . 5 

4.  Results  6 

4.1  V erification  of  the  Model  Insertion . 6 

4.2  Validation  of  the  Model . 7 

4.3  Scalability  of  the  Results . 11 

4.4  An  Application  of  the  Model . 13 

5.  Future  Work  and  Conclusions  25 

6.  References  26 

Appendix.  Microstructural  Equations  28 

Distribution  List  35 


iii 


List  of  Figures 


Figure  1.  Finite-element  mesh  of  a  cylindrical  wedge  under  Taylor  impact . 7 

Figure  2.  DYNA3D  and  AEH-PARDYN  axial  velocity  solutions  for  the  finite-element  mesh 
of  figure  1 . 8 

Figure  3.  DYNA3D  and  AEH-PARDYN  axial  stress  solutions  for  the  finite-element  mesh  of 
figure  1 . 9 

Figure  4.  DYNA3D  and  AEH-PARDYN  effective  plastic  strain  solutions  for  the  finite- 

element  mesh  of  figure  1 . 10 

Figure  5.  Finite-element  model  and  deformed  mesh  used  for  the  AEH-PARDYN  solution  for 
the  verification  problem . 12 

Figure  6.  Finite-element  model  and  deformed  mesh  used  for  the  DYNA3D  solution  for  the 

verification  problem . 13 

Figure  7.  DYNA3D  and  AEH-PARDYN  axial  velocity  solutions  for  the  validation  problem. ...  14 

Figure  8.  DYNA3D  and  AEH-PARDYN  axial  stress  solutions  for  the  validation  problem . 15 

Figure  9.  DYNA3D  and  AEH-PARDYN  effective  stress  solutions  for  the  validation 

problem . 16 

Figure  10.  DYNA3D  and  AEH-PARDYN  effective  plastic  strain  solutions  for  the  validation 
problem . 17 

Figure  1 1 .  Scalability  results  for  a  fixed  mesh  and  for  a  constant  processor  load . 1 8 

Figure  12.  Finite-element  mesh  for  PC74AL37  laminate  and  its  buffer  layer,  back  plate,  and 
flyer . 19 

Figure  13.  Finite-element  mesh  for  PC74SS37  laminate  and  its  buffer  layer,  back  plate,  and 
flyer . 19 

Figure  14.  Finite-element  mesh  for  PC74GS55  laminate  and  its  buffer  layer,  back  plate,  and 
flyer . 20 

Figure  15.  Experimental  result  and  DYNA3D  prediction  for  axial  velocity  at  the  back  of 

buffer  layer  for  the  PC74AL37  laminate . 21 

Figure  16.  Experimental  result  and  DYNA3D  prediction  for  axial  velocity  at  the  back  of 

buffer  layer  for  the  PC74SS37  laminate . 22 

Figure  17.  Experimental  result  and  DYNA3D  prediction  for  axial  velocity  at  the  back  of 

buffer  layer  for  the  PC74GS55  laminate . 22 

Figure  18.  Experimental  result  and  AEH-PARDYN  prediction  for  axial  velocity  at  the  back 
of  buffer  layer  for  the  PC74AL37  laminate . 23 

Figure  19.  Experimental  result  and  AEH-PARDYN  prediction  for  axial  velocity  at  the  back 
of  buffer  layer  for  the  PC74SS37  laminate . 24 

Figure  20.  Experimental  result  and  AEH-PARDYN  prediction  for  axial  velocity  at  the  back 
of  buffer  layer  for  the  PC74GS55  laminate . 24 

Figure  A-l .  Macro/microanalyses  neighborhoods . 28 


List  of  Tables 


Table  1.  Material  properties . 6 

Table  2.  Material  properties . 7 

Table  3.  Wall  clock  times  for  scalability  study . 18 

Table  4.  Laminate  and  impact  velocity  details . 20 

Table  5.  Material  properties . 21 


v 


Intentionally  left  blank. 


vi 


1.  Introduction 


Heterogeneous  materials,  such  as  layered  or  dispersed  particulate  composites,  are  widely  used  in 
military  annor  applications  for  protection  from  ballistic  impacts  or  for  defeating  enemy 
structures  and  vehicles.  For  effective  use  of  these  materials  in  short  duration  impact 
applications,  one  should  be  able  to  quantify  the  local-  or  microlevel  structural  response  of  the 
constituent  phases,  understand  the  different  failure  modes  at  the  microstructural  level,  and 
translate  this  infonnation  into  the  effective  constitutive  material  properties  at  the  global  level. 

Spanning  a  heterogeneous  material  from  its  component  microphases  to  its  global  structural  level 
with  finite  elements  to  obtain  micro-  and  global  structural  responses  leads  to  large  finite-element 
meshes  and  computations  that  are  both  expensive  and  intractable.  Multiscale  analyses  are  used 
to  alleviate  this  problem  and  are  even  being  extended  to  include  nanoscale  analyses  to  aid 
nanomaterial  development. 

At  a  basic  level,  these  analyses  use  one  or  more  microstructural  models  for  describing 
microlevel  material  features,  such  as  fibers,  particles,  lamina,  voids,  etc.,  and  express  the 
microlevel  response  as  a  local  perturbation  of  global  response.  Each  element,  or  a  selected  few 
of  the  global  elements,  calls  the  microstructural  model  to  sample  the  microlevel  structural 
response  and  build  the  element  level  stresses. 

While  many  theories  are  available  for  linking  the  microstructural  response  to  the  global  response, 
the  theory  selected  here  is  capable  of  modeling  elastic/plastic  large  defonnation  seen  in  the 
annor  impact  applications.  It  is  based  on  the  mathematically  rigorous  asymptotic  expansion 
homogenization  (AEH)  method.  While  a  historical  context  to  it  is  presented  in  the  next  section, 
its  most  recent  proponents  are  Chung  et  al.  (1999),  who  extended  its  range  of  application  to 
include  short  transient  loading  and  applied  it  for  predicting  the  elastic/plastic  large  deformation 
of  composite  materials. 

In  the  AEH  approach,  the  microlevel  response  is  expressed  in  terms  of  the  global  response  in  a 
strict,  mathematically  seamless  approach.  An  updated  Lagrangian  scheme  for  small  strains  and 
small  rotations  is  employed  to  account  for  large  displacements,  strains,  and  rotations  over  many 
time  steps.  At  each  time  step,  the  equilibrium  solution  from  the  previous  time  step  is  used  to 
update  the  local  microstructure. 

Since  explicit  time  integrations  take  quite  a  large  number  of  time  steps  and  microlevel 
computations  are  to  be  perfonned  at  each  and  every  time  step  and  for  the  selected  global  finite 
elements,  the  attendant  micro  computations  can  take  much  longer  times  for  solution.  Scalable 
computational  approaches  are  needed  to  use  these  computationally  intensive  methods.  With  the 
aid  of  domain  decomposition  and  message-passing  interface  (MPI),  the  basic  equations  of 
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Chung  et  al.  (1999)  were  recast  for  parallel  computation  and  were  shown  to  be  scalable  for  a 
Taylor  impact  of  a  laminate  with  several  hard  and  soft  layers  by  Valisetty  (2000). 

While  the  previous  two  references  used  the  microstructural  AEH  as  a  stand-alone  code,  the 
current  work  couples  it  as  a  material  level  response  provider  to  PARADYN,  a  general  purpose 
explicit  dynamics  code.  The  objective  was  to  extend  the  applicability  of  AEH  to  large-scale 
multi-material  impact  analyses.  PARADYN  is  a  parallel  version  of  the  Lawrence  Livermore 
National  Laboratory’s  serial  DYNA3D  code.  It  is  a  natural  choice  for  the  present  purpose 
because  of  its  large  number  of  material  models,  contact  algorithms,  element  types,  etc. 

Coupling  two  stand-alone  structural  codes,  such  as  the  microstructural  AEH  and  PARADYN,  is 
best  accomplished  when  the  interface  computer  code  is  written  in  a  form  that  is  as  simple  as 
possible.  While  an  independent  stand-alone  interface  is  desirable,  it  has  to  accomplish  two 
objectives:  (1)  the  data  exchange  between  the  two  codes  must  be  transparent  to  users,  and  (2)  the 
coupling  to  the  driver  code,  i.e.,  PARADYN,  must  be  done  in  such  a  way  that  it  does  not  disturb 
any  of  the  driver  code’s  special  features,  such  as  scalable  parallel  operations,  time  integration, 
contact  algorithms,  etc. 

Since  PARADYN  is  a  scalable  parallel  code,  it  already  has  its  computational  logic  and  data 
structures  set  for  tackling  parallel  computations.  Even  though  the  code  has  facilities  to  add 
user-defined  materials  and  elements,  the  present  microstructural  AEH  is  neither  a  new  material 
nor  an  element.  Since  potentially  large  microstructural  element  models  can  be  considered  with 
the  microstructural  AEH,  requiring  large  amounts  of  microelement  specific  material  properties, 
stresses  and  strains  data  to  be  stored  and  made  available  for  conducting  time  integration,  using 
PARADYN ’s  data  arrays  for  storing  and  retrieving  such  microelement  data  becomes  intractable 
and  potentially  conflicting  to  its  scalability  operations.  To  alleviate  this  situation,  external 
processor  specific  files  are  used  for  writing  a  global  element’s  data  and  microstructural  elements’ 
data.  This  enables  the  AEH  model  coupling  without  disturbing  the  PARADYN  code’s  MPI  and 
domain  decomposition-based  scalability  operations. 

To  enable  access  to  the  microstructure  AEH  model,  a  reference  is  provided  to  it  as  a  “Material 
Type  3”  sub-option  in  PARADYN.  This  reference  is  completely  arbitrary  and  could  be  placed 
even  without  getting  into  any  material  level  subroutines.  Global  elements  tabbed  with  this 
material  option  then  write  the  current  time  step  and  strain  increment  data  to  external  files  and 
read  back  the  effective  mechanical  properties  and  effective  stresses.  With  the  aid  of  a 
user-specified  microelement  finite-element  model  and  after  reading  global  elements’  data,  the 
AEH  model  conducts  microlevel  analyses  and  writes  the  effective  properties  and  stress  data  to 
files.  The  data  exchange  between  the  PARADYN  and  AEH  codes  occurs  through  external 
processor  specific  files.  Thus,  with  this  interface,  a  dual  level  finite-element  modeling  capability 
is  added  to  the  PARADYN  code. 

With  regards  to  writing  the  computer  code  of  the  interface,  while  many  environments  are 
available,  the  one  that  was  considered  is  the  “interdisciplinary  computing  environment  (ICE).” 
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As  described  by  Clarke  and  Namburu  (2002),  this  environment  uses  an  “extensible  data  model 
and  format  (XDMF)”  for  describing  the  nature  of  data  being  exchanged  between  the  codes  and  a 
visual  environment  for  the  user  to  monitor  the  exchanging  data  while  the  codes  run 
independently  in  the  background.  The  xdmf  uses  extensible  markup  language  (XML), 
hierarchical  data  format  version  5  (HDF5),  and  network  distributed  global  memory  (NDGM). 
The  XML  language  uses  PYTHON  language  script  to  enable  users  to  describe  the  types  of  data 
being  exchanged  between  the  independently  running  codes.  The  HDF5  data  format  in  which  the 
codes  exchange  data  via  external  files  is  an  effective  fonnat  not  only  for  doing  read/writes  but 
also  for  grouping  the  data  for  fast  scientific  visualization.  The  NDGM  provides  data  access  to 
codes  running  on  processor  grids  in  a  scalable  manner. 

While  the  ICE  experience  of  data  markup,  visualization,  storage,  and  handling  will  be  discussed 
in  a  future  publication,  the  usefulness  of  the  AEH-PARADYN  coupling  is  presented  in  the 
current  report.  When  coupled  with  PARADYN,  the  microstructural  AEH  can  potentially  avail 
the  various  material  models  for  its  constituent  materials  from  PARADYN.  Also,  in  application 
runs,  PARADYN’s  features  (e.g.,  contact  and  the  ability  to  remove  excessively  deformed  global 
materials)  can  also  be  used.  The  microstructural  AEH  enhancement  of  PARADYN  was 
demonstrated  through  simulation  of  plate  impact  tests  for  three  composite  laminates  made  of 
extremely  dissimilar  materials. 

The  outline  of  this  report  is  as  follows.  Section  2  presents  a  brief  historical  context  to  the  AEH 
development.  Section  3  presents  the  details  of  the  microstructural  AEH  coupling  to  PARADYN 
and  guidelines  for  its  use.  Section  4  presents  verification,  scalability,  and  application  of  the 
coupled  codes.  The  micro-  and  global  structural  response  coupling  and  equations  are  presented 
in  the  appendix. 


2.  AEH  Development 


While  the  rule  of  mixtures  represents  an  earliest  concept  in  composite  micromechanics,  the  idea 
that  stress/strain  distributions  caused  by  holes  and  inclusions  in  an  otherwise  uniform  medium 
can  be  used  to  represent  the  stress/strain  distributions  around  fibers  in  laminas  is  the  starting 
point  for  the  development  of  homogenization  theories.  Licht  Frottement  (1987),  Lene  (1984), 
Moulinec  and  Suquet  (1994),  Bensoussan  et  al.  (1978),  and  Sanchez-Palencia  (1980)  have 
formalized  the  mathematics  of  the  homogenization  theories  for  obtaining  the  effective  material 
properties.  The  works  of  Hashin  (1962),  Hashin  and  Strickman  (1963),  Hill  (1964),  Willis 
(1982),  and  Suquet  (1982)  can  be  seen  as  the  engineering  counterparts  of  these  theories. 

For  a  short  transient  impact  loading  condition,  a  detailed  knowledge  of  the  material  flow  in  and 
around  material  microstructural  constituents  and  material  interfaces  is  necessary.  This  flow  is 
three-dimensional  and  elastic/plastic.  Classical  rule  of  mixtures  theories  based  on  constant  stress 
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or  strain  assumptions  lead  to  an  aggregation  of  the  response  in  the  microstructural  details  and 
thus  are  of  limited  use. 


Recent  advances  in  mathematically  rigorous  AEH  methods  enable  the  coupling  of  micro-  and 
global  approaches  for  both  linear  and  nonlinear  structural  applications  (Fish  et  ah,  1997;  Terada 
and  Kikuchi,  1995;  Ghosh  and  Moorthy,  1995;  Lene,  1986;  and  Guedes  and  Nikuchi,  1990). 
Most  recently,  Chung  et  al.  (1999)  demonstrated  the  applicability  of  the  AEH  method  for 
heterogeneous  media  subjected  to  short  transient  loading  by  employing  explicit  dynamics  finite- 
element  formulations  in  conjunction  with  elastic/plastic  material  response.  The  investigation  of 
Chung  et  al.  (1999)  considers  the  dynamic  equation  of  motion  discretized  according  to  a  second 
order  accurate  scheme  presented  in  Namburu  (1991). 

Using  the  AEH  approach,  the  microlevel  variables  are  expressed  as  direct  functions  of  the  global 
variables  in  a  strict,  mathematically  seamless  approach.  This  method  is  especially  suitable  for 
conducting  explicit  transient  elastic/plastic  analysis  of  heterogeneous  materials.  An  updated 
Lagrangian  scheme  for  small  strains  and  small  rotations  is  employed  to  account  for  large 
displacements,  strains,  and  rotations  over  many  time  steps. 

Because  an  explicit  time  integration  scheme  takes  quite  a  large  number  of  time  steps  and  the 
global  and  microlevel  computations  are  to  be  performed  at  each  time  step  and  for  each  global 
finite  element,  the  attendant  micro  and  global  computations  can  take  much  longer  times  for 
solution.  Scalable  computational  approaches  are  needed  to  explore  these  computationally 
intensive  methods.  In  Chung  et  al.  (1999),  the  microstructural  AEH  was  made  parallel  in  a 
stand-alone  explicit  dynamics  code  with  results  showing  good  scalability  for  the  AEH’s  coupled 
micro/macro  equations.  These  equations  are  presented  in  the  appendix.  They  form  the  basis  of 
the  microstructural  AEH  that  is  now  coupled  with  PARADYN. 


3.  Coupling  the  Microstructural  AEH  Model  in  PARADYN 


PARADYN  is  a  general-purpose  explicit  dynamics  code  with  a  capability  for  defining  multiple 
element  types  and  material  types.  Each  element  of  the  global  finite-element  model  that  is  tabbed 
for  detailed  microlevel  analysis  is  made  to  call  the  microstructural  AEH  model  to  sample  the 
local  microlevel  behavior  and  build  the  element  global  stresses  from  the  constituent 
microstresses. 

3,1  Guidelines  for  Using  the  Model 

The  microstructural  AEH  sub-option  can  be  invoked  by  specifying  a  “Material  Type  3”  and  a 
non-zero  number  in  the  fifth  data  location  on  the  eighth  material  card.  The  material  option 
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describes  a  nonlinear  elastic/plastic  response  with  kinematic  and/or  isotropic  hardening.  The 
added  microcode  is  not  limited  to  microphase  materials  of  this  type,  however. 

The  microstructure  is  next  defined  with  the  usual  nodes,  elements,  and  material  properties. 
Although  any  type  of  format  can  be  used  for  this  purpose,  the  one  that  is  selected  uses  the  same 
PARADYN  input  fonnat.  This  was  done  to  avail  PARADYN’s  preprocessors,  such  as  INGRID. 
Although  not  verified,  at  present  there  are  no  restrictions  on  the  number  of  microstructural 
elements,  nodes,  and  materials. 

All  the  information  defining  the  microstructural  AEH  is  to  be  contained  in  an  input  file  separate 
from  the  PARADYN’s  global  structure-defining  input  file  and  is  to  be  made  available  to  the 
PARADYN  executable  along  with  a  domain  decomposition  part  file. 

3.2  Details  of  the  Data  Exchange 

While  PARADYN  focuses  on  global  computations  such  as  the  aggregation  of  nodal  forces, 
solution  of  the  equations  of  motion,  updating  global  nodal  positions,  enforcing  contact,  etc.,  the 
microstructural  AEH  code  focuses  on  the  microstructural  AEH  calculations  which  typically 
involve  evaluating  the  microstructure  as  a  perturbation  of  global  deformation  and  periodicity. 

For  the  two  codes  to  perform  their  respective  computations,  it  is  necessary  that  the  relevant  data 
be  exchanged  between  them  for  each  time  step  and  for  each  global  element  identified  for 
microcomputations.  In  the  present  work,  the  data  exchange  involves  two  parts  as  follows: 

1 .  From  PARADYN,  at  the  beginning  of  each  time  step  and  for  the  identified  global  elements, 
the  following  data  is  written  to  external  processor-specific  data  files  for  read  access  to  the 
microstructural  AEH  model: 

Current  time  step, 

current  global  strain  increments, 

global  strains  from  previous  time  step, 

global  stresses  from  previous  time  step,  and 

effective  plastic  strain  from  previous  time  step. 

2.  After  reading  the  PARADYN- written  data  files,  the  microstructural  AEH  code  computes 
and  writes  the  following  data  to  external  processor-specific  data  files  for  read  access  to 
PARADYN: 

Effective  density, 

effective  bulk  modulus, 

effective  shear  modulus, 
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computed  current  global  stresses  for  the  element,  and 

computed  current  effective  plastic  strain  for  the  element. 

While  PARADYN’s  native  arrays  are  used  for  tracking  global  stresses,  global  strains,  effective 
properties  (density,  bulk  modulus,  and  shear  strain),  and  effective  plastic  strain,  no  such  native 
arrays  are  used  to  store  the  microelement  quantities.  Instead,  the  microelement  data  are  written 
to,  and  read  from,  processor-specific  files  during  each  time  step  and  updated.  As  opposed  to  the 
global  stresses  of  the  global  elements,  the  stored  values  of  these  responses  reflect  the  accurate 
stress  distribution  within  the  microelements  of  the  microstructure.  During  a  PARADYN  run,  or 
after  the  run  is  complete,  the  processor-specific  files  can  be  visualized  independently  to  analyze 
the  microstructural  stress  and  defonnation  development  for  any  of  the  global  elements. 


4.  Results 


Verification  of  the  coupling  algorithm,  evaluation  of  the  coupling  scalability,  and  validation  of 
the  coupled  codes  were  achieved  through  the  following  numerical  analyses  conducted  on  the 
U.S.  Army  Research  Laboratory’s  SGI  02k  machines. 

4.1  Verification  of  the  Model  Insertion 

For  verifying  the  model  coupling,  a  small  cylinder  with  homogeneous  material  was  considered  in 
a  Taylor  impact,  and  was  analyzed  alternatively  as  a  heterogeneous  material  with  DYNA3D  and 
as  a  homogeneous  material  with  alternating  layers  of  two  different  materials  but  with  same 
identical  properties  using  the  microstructural  AEH  in  PARADYN.  The  idea  was  that  the  AEH 
model  should  default  into  a  homogeneous  material  solution,  and  the  data  exchange  can  point  to 
any  errors  in  the  coupling. 

Enforcing  quarter  symmetry  on  sides  parallel  to  the  axis  of  the  cylinder,  a  90°  wedge  of  the 
cylinder  was  considered  with  108  elements  in  five  rows,  as  shown  in  figure  1.  The  length  and 
radius  of  the  cylinder  was  2.592  mm  and  3.2  mm,  respectively.  The  material  was  considered 
elastic/plastic  with  isotropic  hardening.  The  material  properties  are  shown  in  table  1 . 


Tablet.  Material  properties. 


Density 

(kg/m3) 

Young’s  Modulus 
(GPa) 

Poisson’s  Ratio 

Yield  Stress 
(GPa) 

Tangent  Modulus 
(GPa) 

8930 

117 

0.35 

0.4 

0.1 

The  velocity  of  the  cylinder  was  227  m/s.  To  affect  the  deformation,  one  end  of  the  cylinder  was 
assumed  to  be  free  of  velocity.  Results  for  axial  velocity,  axial  stress,  and  effective  plastic  strain 
are  considered  for  three  nodal  locations,  near  the  free  end,  center,  and  at  the  impacted  end,  as 
shown  in  figure  1 .  They  are  presented  in  figures  2-4,  respectively,  for  the  heterogeneous 
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Figure  1.  Finite-element  mesh  of  a  cylindrical 
wedge  under  Taylor  impact. 

DYNA3D  solution  and  for  the  homogeneous  microstructural  AEH-PARADYN  solution.  Both 
solutions  appear  to  be  in  agreement  with  each  other,  thus  verifying  the  model  coupling  for 
micro  structures  with  homogeneous  materials. 

4.2  Validation  of  the  Model 

In  this  example,  a  cylinder  similar  to  the  one  in  the  previous  example  was  considered  but  with 
two  materials  alternating  in  a  total  of  96  layers  normal  to  the  cylinder  axis.  The  length  and 
radius  of  the  cylinder  were  31.1  mm  and  3.2  mm,  respectively.  The  material  was  considered 
elastic/plastic  with  isotropic  hardening.  The  material  properties  are  shown  in  table  2. 

Table  2.  Material  properties. 


Layer 

Density 

(kg/m3) 

Young’s  Modulus 
(GPa) 

Poisson’s  Ratio 

Yield  Stress 
(GPa) 

Tangent  Modulus 
(GPa) 

1 

8930 

1 17e9 

0.35 

0.4e9 

0.1e9 

2 

8930 

90e9 

0.35 

0.4e9 

0.1e9 

Two  solutions  were  obtained.  The  first  solution  was  a  heterogeneous  DYNA3D  solution 
obtained  using  10,368  global  elements.  To  model  the  repeating  48  sets  of  the  two  dissimilar 
materials,  eight  rows  of  elements  were  used  for  each  material  layer.  The  second  solution  was  the 
homogeneous  microstructural  AEH-PARADYN  solution  obtained  with  a  1296  element  global 
model.  The  elements  are  stacked  in  96  layers.  In  contrast  to  the  first  solution,  a  single  layer  of 
elements  was  used  for  each  material.  All  layers  are  assumed  to  be  exhibiting  homogeneous 
composite  behavior  in  the  global  sense,  but  with  microstructural  AEH  computed  using  a 
two-layer  microstructural  model. 
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Figure  2.  DYNA3D  and  AEH-PARDYN  axial  velocity  solutions  for  the  finite-element 
mesh  of  figure  1 . 
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Figure  3.  DYNA3D  and  AEH-PARDYN  axial  stress  solutions  for  the  finite-element 
mesh  of  figure  1 . 
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Taylor  Test 


Figure  4.  DYNA3D  and  AEH-PARDYN  effective  plastic  strain  solutions  for  the  finite- 
element  mesh  of  figure  1. 
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Enforcing  quarter  symmetry  on  sides  parallel  to  the  axis  of  the  cylinder,  a  90°  wedge  of  the 
cylinder  was  considered.  The  velocity  imposed  on  all  but  one  end  face  of  the  cylinder  is 
227  m/s.  Figures  5  and  6  show  the  two  models  both  before  and  after  the  impact.  Solutions  at 
nodes  near  top,  center,  and  bottom  of  the  cylinder  are  presented  in  figures  7-10,  for  axial 
velocity,  axial  stress,  effective  stress,  and  effective  plastic  strain.  For  all  these  responses,  the 
microstructural  AEH  was  able  to  track  the  DYNA3D’s  heterogeneous  responses  well. 

Deviations  in  the  axial  velocity  and  axial  stress  can  be  attributed  to  three  factors:  (1)  these 
differences  are  near  stress-free  edges  prone  to  be  affected  by  the  transients;  (2)  the  meshes  are 
different  in  their  coarseness  which  affects  the  time  integration  in  the  two  solutions;  and  (3)  the 
values  are  node  averaged  but  the  nodes  are  at  different  locations  in  the  two  models  since  the 
meshes  are  different.  Similarly,  the  differences  in  the  plastic  strain  predictions  can  be  attributed 
to  the  extreme  length  of  the  time  duration  of  the  analysis,  which  makes  the  differences 
accumulate. 

4.3  Scalability  of  the  Results 

During  the  development  of  the  AEH-PARADYN  coupling  interface,  due  consideration  was 
given  to  the  fact  that  PARDYN  is  a  parallel  code  while  coupling  the  microstructural  AEH. 
Because  the  PARADYN  code  has  a  well  tested  MPI/OPENMP -based  computational  logic,  the 
coupling  was  done  by  keeping  a  global  element’s  microstructural  AEH  computations  local  to  the 
processor  on  which  the  global  element  resides.  This  means  that  information  such  as 
microelement  stresses,  strains,  and  mechanical  properties  needed  in  microelement  computations 
are  kept  local  and  not  passed  into  global  data  structures.  Consequently,  a  global  element’s 
microstress,  microstrain,  and  micromechanical  properties  are  needed  to  be  written  to,  and  read 
from,  processor-specific  files  in  between  time  steps.  To  demonstrate  the  fact  that  the  penalty 
associated  with  this  approach  is  small,  results  from  a  limited  scalability  study  are  presented  in 
this  section. 

The  Taylor  impact  model  of  section  4.2  is  used  for  this  purpose.  Two  types  of  scalability 
studies,  one  with  a  fixed  mesh  for  standard  scalability,  and  the  other  with  varying  mesh  size,  but 
with  a  fixed  number  of  elements  per  processor,  are  performed.  The  nodes,  elements,  and  wall 
clock  times  from  these  studies  are  presented  in  table  3.  In  these  simulations,  the  increase  in 
number  of  elements  was  achieved  by  adding  element  rows  in  the  cylindrical  axis  direction  only, 
and  not  in  the  hoop  and  radial  directions. 

The  wall  clock  times  are  also  plotted  in  figure  1 1 .  The  plot  on  the  left  side  of  this  figure  shows 
run  times  for  the  122,880  element  mesh  for  the  standard  scalability  study.  A  linear  speedup  can 
be  discerned  from  this  plot. 

On  the  other  hand,  the  plot  on  the  right  side  of  the  figure  shows  the  run  times  for  the  scaled 
scalability  study.  In  this  study,  the  mesh  size  was  doubled  with  each  doubling  of  the  processors. 
The  run  times  did  not  remain  constant  as  were  the  number  of  elements  per  the  processor.  The 
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Figure  5.  Finite-element  model  and  deformed  mesh  used  for  the  AEF1-PARDYN  solution  for  the 
verification  problem. 

reason  could  be  that  the  plastic  deformation  is  very  concentrated  near  base  and  the  domain 
decomposition  used  to  distribute  the  elements  over  processors  is  not  sensitive  to  this  fact. 
Another  factor  could  be  the  small  processor  load  (-2000  elements/processor).  The  scalability 
might  improve  with  a  higher  processor  load  (e.g.,  20,000  elements/processor  or  higher). 
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Figure  6.  Finite-element  model  and  deformed  mesh  used  for  the  DYNA3D  solution  for  the  verification 
problem. 

4.4  An  Application  of  the  Model 

In  this  section,  the  ability  of  the  microstructural  model  AEH-PARADYN  to  predict  the 
experimentally  observed  impact  behavior  of  laminates  with  widely  dissimilar  and  alternating 
layers  of  materials  is  demonstrated.  For  this  purpose,  impacts  of  three  composite  laminates  from 
the  experimental  work  of  Zhuang  (2002)  are  considered.  In  all  these  impacts,  the  objective  was 
to  study  the  effect  of  composite  lay-ups  with  alternating  hard  and  soft  layers  on  shock  wave 
propagation.  Capturing  the  local  mechanical  behavior  of  such  lay-ups  is  very  difficult, 
especially  in  the  regime  of  shock. 
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Figure  7.  DYNA3D  and  AEH-PARDYN  axial  velocity  solutions  for  the  validation 
problem. 
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Figure  8.  DYNA3D  and  AEH-PARDYN  axial  stress  solutions  for  the  validation  problem. 
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Figure  9.  DYNA3D  and  AEH-PARDYN  effective  stress  solutions  for  the  validation 
problem. 
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Figure  10.  DYNA3D  and  AEH-PARDYN  effective  plastic  strain  solutions  for  the 
validation  problem 
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Table  3.  Wall  clock  times  for  scalability  study. 


Nodes 

Elements 

No.  of  Processors 

Wall  Clock  Time 

9787 

8064 

4 

3072 

18758 

16128 

8 

4051 

35321 

31240 

16 

6046 

71217 

64512 

32 

7492 

137538 

127050 

64 

8944 

140727 

122880 

8 

88369 

140727 

122880 

16 

44376 

140727 

122880 

32 

22303 

140727 

122880 

64 

11353 

Figure  11.  Scalability  results  for  a  fixed  mesh  and  for  a  constant  processor  load. 


In  these  laminates,  a  basic  unit  consisting  of  a  soft  layer  and  a  hard  layer  is  repeated  a  number  of 
times.  Laminates  are  then  bonded  to  soft  buffer  layers  and  mounted  on  to  back  plates.  The  soft 
layer  is  a  polycarbonate  (PC)  plastic  sheet  of  0.74-mm  thickness.  Different  materials  and 
thicknesses  are  used  for  the  hard  layers  as  follows:  aluminum  (0.37-mm  thickness),  steel 
(0.37-mm  thickness),  and  glass  (0.55-mm  thickness).  The  soft  buffer  layer  has  a  thickness  of 
0.74  mm  and  is  made  from  PC.  The  back  plate  has  a  thickness  of  12.5  mm  and  is  made  from 
polymethyl  methacrylate  plastic.  As  shown  in  figures  12-14,  the  laminates  are  impacted  in  a 
direction  normal  to  the  lay-up  using  flyer  plates  made  of  PC.  The  flyer  plate  has  a  thickness  of 
2.87  mm.  The  diameters  of  the  flyer  plates,  laminates,  and  the  back  plates  are  34,  38.1,  and 
38.1-mm,  respectively.  The  details  of  the  three  laminate  constructions  are  in  table  4. 
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Figure  12.  Finite-element  mesh  for  PC74AL37 

laminate  and  its  buffer  layer,  back  plate, 
and  flyer. 


Figure  13.  Finite-element  mesh  for  PC74SS37 
laminate  and  its  buffer  layer,  back 
plate,  and  flyer. 
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Figure  14.  Finite-element  mesh  for  PC74GS55 
laminate  and  its  buffer  layer,  back 
plate,  and  flyer. 


Table  4.  Laminate  and  impact  velocity  details. 


Laminate  Types 

No.  of  Repeating 
Units 

Laminate  Thickness 
Including  Buffer  Layer 
(mm) 

Flyer  Plate  Velocity 
(m/s) 

PC74AL37 

5 

6.62 

589 

PC74SS37 

8 

9.97 

561 

PC74GS55 

7 

9.95 

563 

As  can  be  seen  from  the  deformed  meshes  in  figures  12-14,  different  meshes  are  used  for  the 
three  laminates.  The  PC74AL37  laminate,  including  its  buffer  layer,  back  plate,  and  the  flyer 
plate,  were  all  meshed  with  a  total  of  69899  nodes  and  55095  elements.  The  PC74SS37  laminate 
and  buffer,  back,  and  flyer  plates  were  meshed  with  a  total  of  108326  nodes  and  82827  elements. 
Finally,  the  PC74GS55  laminate’s  impact  was  modeled  with  a  total  of  167549  nodes  and  141798 
elements.  Because  the  experiments  showed  much  delamination,  and  because  the  interlaminar 
strengths  of  the  laminates  are  relatively  small,  i.e.,  when  compared  to  the  impact  generated  loads, 
laminar  interfaces  were  used  and  modeled  as  zero  strength  slide  surfaces.  The  effect  of  these 
slide  surfaces  can  be  seen  as  outflowing  layers  in  the  defonned  meshes  in  figures  12-14.  The 
material  properties  used  for  the  present  simulations  are  in  table  5. 
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Table  5.  Material  properties. 


Material 

Density 

(g/m3) 

Shear  Modulus 
(GPa) 

Yield  Stress 
(GPa) 

Plastic  Modulus 
(GPa) 

Poisson’s  Ratio 

PC 

1190 

0.94 

0 

1.60 

0.37 

PMMA 

1180 

1.20 

0 

1.60 

0.34 

6061  AL 

2710 

30.0 

0.32 

0.69 

0.33 

304  SS 

7890 

77.0 

0.33 

1.7 

0.29 

D  263  glass 

2510 

30.1 

— 

— 

0.208 

Before  conducting  the  dual  level  macro/micro  analyses  with  AEH-PARADYN,  single  level 
DYNA3D  global  analyses  were  conducted  using  the  aforementioned  finite-element  models. 
Using  DYNA3D’s  elastic/plastic  material  Model  3,  for  layer  materials,  the  laminate  is  treated  as 
materially  heterogeneous.  The  results  for  the  axial  velocity  at  the  back  of  the  buffer  layer  are 
compared  with  the  experimental  results  in  figures  15-17.  Even  though  the  overall  trends  of 
loading  and  unloading  of  the  shock  wave  are  predicted  well,  there  is  much  discrepancy  in  the 
predictions  of  the  initial  slope  and  the  peak  velocity  reached. 


Figure  15.  Experimental  result  and  DYNA3D 

prediction  for  axial  velocity  at  the  back 
of  buffer  layer  for  the  PC74AL37 
laminate. 
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Figure  16.  Experimental  result  and  DYNA3D 
prediction  for  axial  velocity  at  the 
back  of  buffer  layer  for  the 
PC74SS37  laminate. 


Figure  17.  Experimental  result  and  DYNA3D 
prediction  for  axial  velocity  at  the 
back  of  buffer  layer  for  the 
PC74GS55  laminate. 

The  AEH-PARADYN  analyses  were  run  assuming  that  the  laminates  were  all  homogeneous  in 
the  global  sense,  but  heterogeneous  at  microlevel.  Two-layer  microstructural  AEH  models  were 
used  to  bring  in  the  lamination  effect  of  the  hard  and  soft  layers.  At  the  microlevel,  the  layer 
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constituents  are  assumed  to  exhibit  elastic/plastic  behavior  that  is  similar  to  the  DYNA3D’s 
elastic/plastic  Material  Type  3.  The  predictions  from  these  analyses  for  the  axial  velocity  at  the 
back  of  the  buffer  layer  are  compared  with  the  experimental  results  in  figures  18-20.  For  all  the 
three  laminates,  the  predictions  compared  well,  especially  in  following  the  initial  slope,  peak 
velocity,  and  the  post  peak  unloading. 

This  example  clearly  demonstrates  the  ability  of  the  AEH  microstructural  model  to  accurately 
predict  the  shock  wave  propagation  in  highly  dissimilar  layered  materials,  and  the  capability  of 
the  included  AEFI  microstructural  model  for  conducting  multiscale  macro/micro  analyses  with 
AEH-P  ARAD  YN . 
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Figure  18.  Experimental  result  and  AEH- 
PARDYN  prediction  for  axial 
velocity  at  the  back  of  buffer 
layer  for  the  PC74AL37 
laminate. 
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Figure  19.  Experimental  result  and  AEH- 
PARDYN  prediction  for  axial 
velocity  at  the  back  of  buffer  layer  for 
the  PC74SS37  laminate. 


Figure  20.  Experimental  result  and  AEFI- 
PARDYN  prediction  for  axial 
velocity  at  the  back  of  buffer 
layer  for  the  PC74GS55 
laminate. 
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5.  Future  Work  and  Conclusions 


Microstructures  are  features  of  the  U.S.  Army’s  advanced  structures  and  annor  materials. 
Homogenization  theories  such  as  the  AEH  are  being  developed  to  model  the  effects  of  these 
microstructures  on  global  material  and  structural  responses.  With  the  present  coupling  of  the 
AEH-based  microstructural  model  in  PARADYN,  many  in-situ  material  responses  can  be 
studied  under  diverse  global  structural  conditions.  The  sample  applications  presented  in  this 
report  demonstrate  this.  Efforts  underway  include  the  following  extensions  to  this  work: 

•  using  PARADYN’s  diverse  material  library  to  model  different  microphase  materials, 

•  modifying  the  AEH  model  used  in  this  work  to  include  an  equation  of  state,  and 

•  extending  the  application  range  to  include  complex  microstructures  and  diverse  global 
applications. 
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Appendix.  Microstructural  Equations 


A.l  Governing  Equations 

These  equations  are  presented  in  several  references,  most  recently  in  Chung  et  al.1  Assume  a 
three-dimensional  body  Q  is  an  assembly  of  periodic  structures  containing  different  materials,  as 
shown  in  figure  A-l.  Typically,  the  unit  cell  is  very  small  and  of  order  s  (where  s  is  a  small 
positive  number)  compared  to  the  dimensions  of  the  problem  domain.  Structural  response 
quantities  such  as  displacements,  velocities,  stresses,  and  strains  are  assumed  to  have  slow 
variations  (macroscopic)  from  point  to  point  as  well  as  fast  (microscopic)  variations  within  a 
small  neighborhood  s  of  a  given  point  x.  Let  Y  be  a  (periodic)  representative  part  of  Q. 


Here  we  distinguish  two  scales:  the  macroscopic  scale  (x  e  Q)  and  microscopic  scale 
(yefi).  Lagrangian  representations  of  conservation  of  mass  and  momentum  equations  are  used. 
The  conservation  of  mass  is  used  to  calculate  the  current  density  from  the  initial  density.  The 
momentum  equation,  the  kinematic  relations,  and  the  constitutive  model  are  represented  by 


pV.-a’u-f,, 

(A-l) 

^A(<.+Vy, 

(A-2) 

and 


'Chung,  P.  W.;  Namburu,  R.  R.;  Tamma,  K.  K.  Three-Dimensional  Elasto-Plastic  Heterogeneous  Media  Subjected  to  Short 
Transient  Loads.  AIAA  paper-99-1239  1999. 
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.  8 


<7 ;=Sv(ejtE), 

(A-3) 

where  p  is  the  density,  f.  is  the  body  force  vector,  E  is  the  internal  energy,  the  superscript  s 
denote  micro/macro  continuum  solutions,  v,  is  the  velocity  vector,  cr.  is  the  stress  tensor,  and  e(/ 

is  the  strain  tensor.  The  initial  conditions  in  the  domain  Q '  and  boundary  conditions  on  the 
surface  of  the  domain  f  ( =  STj  u  dT  2 )  are  given  by 

uf(t  =  0)  =  u£Oj, 

(A-4) 

v,'(f  =  0)  =  v', 

(A-5) 

u';  =  if  on 

(A-6) 

and 

ajini  =  Tt  on  Sr2, 

(A-7) 

where  ut  is  the  displacement  vector  and  nt  is  a  surface  normal. 


A.2  Time  Integration 


A  second-order  accurate  Lax-Wendroff-based  explicit  time  integration  procedure  Tamma  and 
Namburu2  is  employed  for  the  conservation  equation  or  equation  of  motion.  In  this  approach, 
the  dependent  variable  velocity  is  first  discretized  in  time  using  a  second  order  Taylor  series 
expansion. 


+  A  tvf 


r+1/ 2 


(A-8) 


Taking  appropriate  time  derivatives  of  equation  1  and  substituting  for  the  velocity  and 
acceleration  gives 


Pvi 


pv, 


=  \<J 


jhj 


+  Pf) 


n+ 1/  2 


(A-9) 


The  stress  increment  is  related  to  a  strain  increment  through  an  appropriate  constitutive  equation. 
The  elasto-plastic  constitutive  equation  at  midpoint  time  increment  is  defined  by 


n— 1/2  1/2 

+  CU  £s  At' 

klmn  mn 


+1/2 


)> 


(A- 10) 


where  Mjjkl  and  Cm  are  radial  return  in  plasticity  tensor  and  elasticity  tensor. 


Stresses  and  strain  rates  at  midpoint  time  can  be  evaluated  using  the  predicted  displacements  as 
shown  by 


2 

Tamma,  K.  K.;  Namburu,  R.  R.  A  Robust  Self  Starting  Explicit  Computational  Methodology  for  Structural  Dynamic 
Applications:  Architecture  and  Representations.  International  Journal  for  Numerical  Methods  in  Engineering  1990,  29,  1441  — 
1454. 
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U; 


r+1/2 


+1/2 


(A-l  1) 


A.3  Asymptotic  Expansion  and  Elasto/plastic  Constitutive  Equations 

The  homogenization  method  is  based  on  the  asymptotic  expansion  of  the  primary  variables 
together  with  a  unit  cell  approach  for  a  heterogeneous  structure.  Assume  a  3-D  body  Q  is  an 
assembly  of  periodic  structures.  Typically,  the  unit  cell  is  very  small,  of  order  <?(£is  a  small 
positive  number)  compared  to  the  dimensions  of  the  problem  domain.  Typically,  two-scale 
asymptotic  expansion  can  be  employed  to  approximate  the  displacement,  velocity,  strain,  or 
stress  fields. 

u ■  (x,  y )  =  uQj  (x, y )  +  sux  (x, y )  +  s2u2  (x,  y)  +  •  •  • ,  (A-l 2) 


vf  (x,  y )  =  v0j  (x,  y)  +  svx  (x,  y)  +  s2v2j  (x,  y )  +  •  •  • , 


(A- 13) 


< (x,  v)  =  e0  (x, y )  +  seh  (x, y)  +  £2e2j  (x,  y )  +  •  •  • ,  (A- 1 4) 

and 

< (A  v)  =  cro  (x, y)  +  £-cr,  (x, y)  +  £2<t2  (x, y)  +  •  •  • ,  (A-l 5) 


where  u‘(x,y) ,  vf(x,y),  eEt(x,y) ,  and  cr/  (x,  y)  are  Y-periodic  functions. 


The  asymptotic  expansion  homogenization  (AEH)  approach  for  nonlinear  applications  is  based 
on  the  instantaneously  linearized  assumption  for  the  constitutive  model.3  The  fundamental 
assumption  of  the  AEH  approach,  therefore,  is  that  the  true  solution  in  the  s  space  is  decomposed 
into  a  macro  space  x  and  a  micro  space  y.  The  basic  assumption  here  is  that  multiple  scales  exist 
only  in  the  spatial  variables  and  that  no  such  scaling  exists  for  the  time  variable.  Further  reading 
in  space  and  time  asymptotic  expansion  approaches  can  be  found  in  Bensoussan.4  The  general 
approach  to  heterogeneous  problems  is  to  separate  and  draw  clear  distinction  between  the  micro- 
and  macrolevel  equilibrium  equations  regardless,  per  the  earlier  linearized  assumption,  of 
material  nonlinearity.  This  is  accomplished  by  asymptotically  expanding  the  primary  variables 
where  the  asymptotic  scale  is  approximated  to  the  second  order. 

vf  (x,  y)  =  vo  (x,  y )  +  svx  (x,  y)  +  s2v2j  (x,  y )  +  •  •  • .  (A- 1 6) 


Spatial  gradients  in  s-space  are  taken  with  respect  to  the  x-coordinate  system.  The 
corresponding  gradient  for  Y-periodic  functions  (where  micro  and  macro  are  now 
distinguishable)  is  given  by  the  chain  rule  where  the  scaling  is  defined  by  y  =  x/sr,  hence,  for  any 
Y-periodic  function  (j),  the  earlier  gradients  are  replaced  by 


Terada,  K.;  Kikuchi,  N.  Nonlinear  Homogenization  Method  for  Practical  Applications,  Computational  Methods  in  Micro¬ 
mechanics.  American  Society  of  Mechanical  Engineers,  Applied  Mechanics  Division,  1995,  212,  1-16. 

4Bensoussan,  A.;  Lions,  J.  L.;  Papanicolaou,  G.  Asymptotic  Analysis  for  Periodic  Structures;  New  York:  North  Holland 
Publishing  Co.,  1978. 
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(A- 17) 


d(f)  _d(j)  1  d(f> 
dx,  fix,  a  dy, 


Then,  using  equation  A- 1  in  equation  A-7  while  considering  equation  A-2  gives  the  rate  of  strain 
tensor  defined  by 

%  =  K  <y0 )  +  4  (V1 )]  +  4^  (fi  )  +  4  ( v2  )]+••• ,  (A- 1 8) 


where  the  symmetrized  gradient  tensors,  ex ,,  and  e’(/  are  defined  by 

1  Sfi  dd>, 

=  -(—  +  —) 

2  fix,  fix. 


and 


4W=i(M+M). 

’  2  8r,  V 


(A- 19) 


(A-20) 


The  relevant  expressions  for  gradients  and  velocities  are  first  substituted  into  the  governing 
equations  of  motion.  Next,  the  micro  and  macro  equations  are  identified  by  selecting  the 
appropriate  coefficients  to  the  scaling  factors  s  which  must  each  be  identically  zero.  Substituting 
the  asymptotically-expanded  velocities,  equation  A-l,  in  the  strain  rate,  equation  A-3,  and  then 
using  the  constitutive  equation,  equation  A- 16,  in  the  equations  of  motion,  equation  A- 15,  yields 
a  set  of  equations  dependent  on  powers  of  e.  To  satisfy  the  equations  of  motion,  each  tenn 
associated  with  each  of  the  powers  of  £  must  approach  zero  identically.  This  leads  to  a  set  of 
equations  associated  with  the  microscopic  and  macroscopic  equations  of  motion.  The  first 
equation,  associated  with  the  powers  of  s  2,  is  given  by 


-2 


dv 


n- 1 


dy,Cm  dy, 


=  0, 


(A-21) 


where  C*jkl  =MijpqCpqU.  Equation  A-21  states  that  v"'1  is  a  function  only  of  x.  Hence, 

derivatives  of  y  are  zero.  Using  this  inference,  the  equation  associated  with  the  powers  of  £  !  is 
derived  as 


C 


dv 


2-1 

ok 


ijkl 


■  +  ■ 


C 


fiv, 


2—1 


i  a 


fix,  tyj  l]kI  dy,  At  dVj 


(Mijklank-V2). 


(A-22) 


Equation  A-22  relates  the  perturbative  velocity  field,  vi,  to  the  macroscopic  velocity  field,  vo. 

The  equation  relating  these  two  quantities  will  become  the  important  micro/macroequation, 
providing  the  direct  link  to  relate  microscopic  to  the  macroscopic  velocities.  It  provides  the 
corrections  which  account  for  the  shape  of  the  interface  separating  various  phases  and  the  time- 
dependent  plastic  softening  effect  due  to  material  nonlinearities.  For  a  homogeneous  material 
problem,  where  gradients  over  y  of  C*/dov"  1  /  fix,  are  zero,  and  in  which  plastic  effects  have  noy 
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dependence,  the  solution  to  equation  A-22  v"  1  is  identically  zero.  In  the  special  case  where  no 
plastic  yielding  occurs,  the  problem  degenerates  to  a  transient  elastic  homogenization  problem. 

Finally,  the  equation  associated  with  the  powers  of  s°  is  written  as 


3v. 

a*.  m 


n- 1 


+— c; 

dy 


dx, 

dv"~l 


+  —C' 


dv‘ 


n- 1 


ijkl 


dx. 


d  * 

+— Cl 

dy 

d 


ijkl 


dy, 

dv"-1 

zk 

dy, 


A,  dx.  Ms"CT" 


P  rn+ 1/2 
A  .  J  i 
At 


(A-23) 


Much  like  equation  A-22  provided  the  basis  for  the  microscale  problem,  equation  A-23  provides 
the  basis  for  the  macroscale  problem.  By  taking  the  volume  average  of  equation  A-23, 
substituting  for  v,  using  the  solution  of  equation  A-22  and  noting  the  periodicity,  the  final 

governing  equation  of  motion  is  given  by 


a 


p)Av:;1  =  A^«1/2  +^)  +  A  t(p)f 


-n+l/2 


(A-24) 


where  the  “corrector  stress,”  a,, ,  is  defined  by 


r~c  _ 

a  ij  ~  ^  ijkl 


dv[ 


n- 1 


dy, 


and  where  quantities  in  brackets  denote 


« )dY ■ 


(A-25) 


(A-26) 


It  is  of  interest  to  note  that  equation  A-24,  for  heterogeneous  conditions,  is  similar  in  form  to 
equation  A-6,  the  homogeneous  equation,  making  the  integration  of  a  micromechanical  problem 
into  a  macroanalysis  tractable  and  straightforward. 

A.4  Constitutive  Equation:  Elasto/Plastic 

The  stress  tensor  is  split  into  two  parts.  The  first  part,  Sy,  is  deviatoric  stress,  which  is  related  to 
material  strength,  and  the  second  part  is  pressure,  P.  To  simplify  notation,  the  superscript  ^has 
been  omitted. 

av=-SyP+Sy  (A-27) 

and 

=~P-  (A-28) 
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The  first  term  in  equation  A-27  accounts  for  volumetric  changes  and  is  typically  evaluated  from 
equations  of  state  in  explicit  dynamics.  The  second  term,  the  stress  deviator,  is  related  to  the 
deformation  of  the  material  and  is  defined  by  a  constitutive  model,  and  in  particular,  the  time  rate 
of  change  of  the  deviator  is  evaluated  from  the  strain  rates. 

5//=26V;,+A,.  (A-29) 

This  model  can  be  defined  with  any  number  of  local  nodes,  elements,  and  materials,  but  since 
this  model  is  called  repeatedly  for  each  time  step,  for  each  element  for  which  microstructural 
response  is  deemed  important,  a  complex  microstructure  model  can  overwhelm  the  overall 
computations,  therefore  keeping  the  model  simple  is  important.  A?/  is  the  correction  for  rigid 
body  rotation. 

An  isotropic  hardening  model  with  a  rate-dependent  Von-Mises  yield  condition  is  employed  in 
the  present  analysis.  A  consistency  condition  ensures  that  the  stress  state  remains  on  the  yield 
surface  at  the  start  and  end  of  a  time-step.  Such  a  condition  is  given  by 

Kn+1Qy=S;+1.  (A-30) 

The  variable  R  is  the  rate-dependent  radius  of  the  yield  surface  or  the  apparent  yield  stress,  Qn  is 
the  vector  specifying  the  normal  direction  to  the  yield  surface,  and  S:j  is  the  deviatoric 
component  of  the  stress.  The  stresses  are  understood  to  be  the  co-rotated  second-order  tensor 
according  to  the  Jaumann  definition  of  the  co-rotational  derivative.  The  Jaumann  derivative  of 
the  Cauchy  stress,  oi  y,  is  related  to  the  material  time  derivative,  a  y,  by 

&ij  =  &y  +  alkajkj  +  ojlkak/ ,  (A-3 1 ) 


where  coy  is  the  rotation  tensor,  the  standard  skew-symmetric  component  of  the  velocity  gradient. 
The  radius  of  the  yield  surface  is  defined  by 


R  =  crl 


1  + 


f  e  ^ 


D) 


+  H'sp , 


(A-32) 


where  d,,  is  the  static  yield  stress  of  the  material,  and  D  and  p  are  the  so-called  fluidity 

—  p 

parameters.  H'  is  the  hardening  parameter,  and  s  is  the  effective  plastic  strain.  The  effective 
rate  of  deformation,  c ,  is  defined  by 


e  = 


(A-3  3) 


where  ey  is  the  deviatoric  component  of  the  total  strain,  Sij.  The  normality  condition,  specified 
by  the  normal  vector  Qy,  is  given  by 


Qij  = 


(A-34) 
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For  a  homogeneous  material  (which  is  applicable  in  the  present  derivation  since  the  constitutive 
equations  are  applied  at  the  microlevel  where  the  material  is  homogeneous  within  each  phase), 
the  incremental  relationships  for  the  apparent  yield  stress  and  the  deviatoric  stress  are  defined  as 


R"+1  =R"+-H'AA 


(A-35) 


and 

S;+l  =  sf  -2GAAQV,  (A-36) 

T 

where  S"+l  is  the  deviatoric  trial  stress,  A/l  is  a  scalar  quantity  representing  the  magnitude  of 

y 

the  radial  return  correction  due  to  plastic  yielding,  and  G  is  the  elastic  shear  modulus.  The 
deviatoric  trial  stress  is  defined  by 

Sf  =S;+C ^  At,  (A-37) 


where  C(/vc/  is  the  tensor  containing  elastic  properties. 


Finally,  substituting  equation  A-37  and  A-35  into  equation  A-30  and  solving  for  the  radial  return 
correction  parameter  gives 


AT  = 


3/2 


c 3G  +  H ) 


n+ 1  o«+ 1 


Sr  -Rn 


(A-38) 


Using  equation  A-38  in  equation  A-36  and  employing  the  standard  assumption  that  a  material  is 
elastic  in  its  dilatational  behavior  and  plastic  only  in  shear  gives  the  total  stress  increment  in 
equation  A- 16.  The  constitutive  equation  A- 16  can  now  be  employed  in  the  derivation  of  the 
micro  and  macro  governing  equations. 
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